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ABSTRACT 

In this study we explore the nucleosynthesis in the dynamic ejecta of compact binary 
mergers. We are particularly interested in the question how sensitive the resulting 
abundance patterns are to the parameters of the merging system. Therefore, we sys- 
tematically investigate combinations of neutron star masses in the range from 1.0 to 
2.0 Mq and, for completeness, we compare the results with those from two simulations 
of a neutron star black hole merger. The ejecta masses vary by a factor of five for the 
studied systems, but all amounts are (within the uncertainties of the merger rates) 
compatible with being a major source of cosmic r-process. The ejecta undergo a ro- 
bust r-process nucleosynthesis which produces all the elements from the second to the 
third peak in close-to-solar ratios. Most strikingly, this r-process is extremely robust, 
all 23 investigated binary systems yield practically identical abundance patterns. This 
is mainly the result of the ejecta being extremely neutron rich (Y e m 0.04) and the 
r-process path meandering along the neutron drip line so that the abundances are 
determined entirely by nuclear rather than by astrophysical properties. While further 
questions related to galactic chemical evolution need to be explored in future stud- 
ies, we consider this robustness together with the ease with which both the second 
and third peak are reproduced as strong indications that compact binary mergers are 
prime candidates for the sources of the observed unique heavy r-process component. 
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1 INTRODUCTION 

About half of the elements heavier than iron are formed 
by neutron capture reactions that occur rapidly in com- 
parison with /3-decays. The basic physical mechanisms of 
this "rapid neutron capture" or "r-process", for short, had 
been identified already in the seminal paper of IBurbidge] 
et al. (19571. Nevertheless, the cosmic cauldrons in which 



these heavy elements are forged have remained elusive for 
more than half a century. Observations of metal-poor stars 
point to at least two groups of r-process events. The first 
one is rare and produces whenever it occurs the heaviest 
r-process elements in nearly exactly solar proportions. The 
second one, occurs more frequently and produces predomi- 
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nantly lighter elements from strontium to silver (Cowan & 
Sncdcn 2006 Hond a~et al.|2006 |. Given that its signature is 



less unique the second component might also be the result 



of a superposition from several sources (see e.g. Arcones & 
|Montes| ( poTT] ) ; |Frischknecht et al] ( |2012[ )) . Metal-poor stars 
that are enriched by r-process material show patterns that 
very closely match the (scaled) solar system abundances for 
nuclei beyond Ba (Z=56). This suggests that this heavy r- 
process component is produced only if a very unique set of 
astrophysical conditions is realized, or, alternatively, that it 
is produced in a range of realizations, but is insensitive to 
the exact parameters of its formation environment. So far, 
there is no generally accepted explanation for the unique 
heavy r-process component. 

Historically supernovae were considered the "natural" 
source of r-process elements, but this view has recently been 
distressed by a slew of investigations (e.g 



p007| |; |Roberts et al.| ( |2010| ); |Fischer et al.| ( |2010| ; |Hudepohl 



Arcones et al. 
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Figure 1. Sensitivity to mass ratio: shown are density cuts ([600 km x 600 km], colour-coded is the logarithm of density in cgs units) 
of a 1.2-1.4 M (t= 13.9 ms), a 1.4-1.4 M© (t= 13.4 ms) and a 2.0-1.4 Mq (t= 15.0 ms) merger. 



et al. (201CM). These studies found neutrino-driven winds 
in supernovae to be seriously challenged in providing the 
physical conditions (high entropy, low electron fraction to- 
gether with rapid expansion) that are required to produce 
the heavy (A > 90) r-process elements. A possible exception 
may be magnetorotationally driven supernova jets where the 
fast ejection of highly compressed neutron-rich matter pro- 



duces interestingly low electron fraction values (Winteler 
|et al.|[2012[ ). It remains to be explored, however, how ro- 
bust this scenario is with respect to the stellar parameters 
and with respect to its nucleosynthetic yields. The main 
contenders of supernovae in terms of r-process nucleosyn- 
thesis are compact binary mergers of either two neutron 
stars (ns 2 ) or a neutron star and a stellar-mass black hole 
(nsbh) (jLattimer fc Schramm|1974l[l976lpc"hler et al.|1989 



Freiburghaus et al. 1999 1 



Since the first nucleosynthesis calculations for neutron star 
merger ejecta (Freiburghaus et al. 19991 have been per- 



formed, numerical models have seen vast improvements in 
physics, numerical methods and resolution and we consider 
it timely to address this important topic again. There have 
been recent r-process studies in a compact binary merger 
context ( Roberts et al.|2011 Goriely et al.|20lT Wanajo & 



Janka|2012 \ which, in agreement with the results of the first 
study, confirm the occurence of a robust r-process. While all 
previous studies provided new insights each of them had its 
limitations. For example, in all of them the ejecta electron 
fraction was poorly known since y e -changing weak interac- 
tions were ignored. More importantly, the ns 2 parameter 
space was explored only very punctually around ns masses 
of 1.4 Mq ( Freiburghaus et al. ( 1999 1 explored one case from 



Rosswog et al. (19991, two cases were explored in the con- 
formal flatness approximati on approach of Goriely et al.| 
( 201l| p| and three cases in 



Roberts et al. 



(20111). In re- 



cent years, however, it has become clear that the neutron 
stars realized in nature span a substantial range in masses, 
see e.g. Lattimer & Prakash (20101 for compilation of neu- 



tron star masses. Therefore, an important open question 
needs to be answered: "How sensitive is the compact object 



1 These authors also used the Lattimer-Swesty EOS to test for 
the impact of the neutron star equation of state. 



merger nucleosynthesis to the actual astrophysical system 
that merges?" 

This is what we address in this paper: we systematically scan 
the neutron star binary parameter space in 21 simulations 
and we calculate the resulting nucleosynthetic products. For 
completeness, we also show the results for the ejecta of two 
nsbh systems. Apart from exploring a very large neutron star 
mass range, we also significantly improve on the treatment 
of the electron fraction Y e : we start with cold neutron stars 
in /3-equilibrium and allow for the evolution of Y e via elec- 
tron and positron captures (Rosswog fc Liebendorfer| 2003). 
We demonstrate that the nucleosynthesis results are essen- 
tially independent of the parameters of merging system, i.e. 
every ns 2 or nsbh system produces nearly exactly the same 
abundance pattern. The results show, however, some sensi- 
tivity to the properties of extremely neutron-rich nuclei (e.g. 
binding energy, half-lives, fission properties) which are not 
well known ( Thielemann et al.|20Tl l. 

Although heavy element nucleosynthesis is the clear focus of 
this study, we would like to stress the importance of compact 
object merger ejecta and their nucleosynthesis in a broader 
astrophysical perspective. Both the LIGO and Virgo gravi- 
tational wave (GW) detectors are currently being upgraded 
( Abbott et al.|2009|[Smith|2009|[Sengupta, et al.|2010| ) to a 
~ 10 — 15 times better sensitivity than the original versions 
of the instruments. This will increase the volume of acces- 
sible astrophysical sources by more than a factor of 1000 
reaching a detection horizon of a few hundred Mpc for ns 2 



mergers and about a Gpc for nsbh mergers (Abadie et al 
2010 1. Since the first detections will most likely be around or 



even below threshold it is of paramount importance to iden- 
tify accompanying electromagnetic signatures of the most 
promising GW sources to enhance the confidence in a possi- 
ble detection. Short gamma-ray bursts are likely the bright- 
est electromagnetic manifestations that result from a com- 
pact binary merger, but also the r-process nucleosynthesis 
and the decay of the resulting radioactive nuclei may leave 
an observable signature as optical/UV transients, so-called 



"macronovae" 


(|Li & Paczyriski||1998 Kulkarni||2005 


Ross- 


wog| 2005 


|Metzger et al.|2010||Roberts et al.|2011||Metzger 


& Berger 


2012 


1. For neutron star mergers they are expected 



to peak about 10 hours after the mergers with a few times 
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10 42 erg/s (Rosswog et al. 



2012 



Piran et al. 



20121. Later, 



on time scales of months to years, the ejecta dissipate their 
kinetic energy in the ambient medium and thereby produce 
possibly detectable radio flares ( |Nakar fc Piran|2011 Piran 
et al.|2012 i. 

This paper is structured as follows. In Sec. [2] we briefly sum- 
marize our simulations and explore the amount of ejected 
matter as a function of the binary system parameters. In 
Sec. [3] we describe how we calculate the nuclear abundances 
and we discuss how they depend on the astrophysical prop- 
erties of the merging system and on the nuclear physics in- 
put. In Sec. [4] we summarize our results and discuss their 
astrophysical implications. 



2 HYDRODYNAMIC EVOLUTION AND 
DYNAMICALLY EJECTED MASS 

2.1 Explored parameter space 

Until recently, the neutron star mass distribution 
thought to be clustered narrowly around 1.35 M (|Thorsett 



& Chakrabarti 1999). Therefore, essentially all neutron star 
merger studies have focused on a narrow range of masses 
around this value. Recent observations, however, indicate a 
much broader neutron star mass spectrum. There is now 
ample support for neutron star masses significantly larger 
than 1.5 M@. A broad peak around 1.5-1.7 Mq has been 



found (Kiziltan et al. 2010 Valentim et al. 20111 for neu- 



tron stars with white dwarf companions, an additional low- 
mass peak near 1.25 Mq is thought to be characteristic for 
neutron stars that were produced by electron capture su- 
pernovae ( Podsiadlowski ct al. 2004 van den Heuvel||2004 



Schwab et al.|2010| >. PSR J1614-2230 with 1.97 ± 0.04 M 

possesses the largest accurately known neutron star mass 
( Demorest et al.|2010 1 , but even higher masses may exist in 
nature. For example, the mass of the black widow pulsar has 
recently been determined as 2.4 ± 0.12 M0(van Kerkwijk 



0.12 M _ 

et al.|20TT l. On the lower mass side, the secondary neutron 
star in J1518+4904 has a best estimated mass value of only 
0.72 Mq, although with a very large la error bar of 0.5 Mq, 
see Lattimer & Prakash (20101 and references therein. 
We take these results as a motivation for a systematic ex- 
ploration of the parameter space from 1.0 to 2.0 Mq in steps 
of 0.2 Mq. Since the neutron star viscosity cannot substan- 
tially spin up the neutron stars during the short tidal in- 



teraction phase preceeding the merger (Bildsten & Cutler 



1992 Kochanek 



1992) all our models have vanishing initial 



ns spins. Our ns cases are complemented by two nsbh cases 
with black hole masses of 5 and 10 Mq, for an overview over 
the performed simulations see Tab. [I] 



2.2 Methodology 

To follow the hydrodynamic evolution we use the Smooth 



Particle Hydrodynamics (SPH) method, see Monaghan 



(2005} and |Rosswog| ( |2009 |) for recent reviews, which, due to 



its completely Lagrangian nature, is ideally suited to follow 
the ejected material. Our code is an updated version of the 



one that was used in earlier studies ( Rosswog & Davies 2002 



Rosswog & Licbendorfcr 2003 , Rosswog ct al. 2003 . Rosswog 
20051. It uses the Shen et al. equation of state (EOS) (Shen 



et al.|1998a|b \ and an opacity-dependent multi-flavour neu- 
trino leakage scheme (Rosswog & Liebendorfcr 2003). The 
latter allows in particular to follow the Y e -changes due to 
electron and positron captures. Moreover, we use a modern, 
time-dependent artificial viscosity prescription, see [Rosswog] 
et al. (2000 20081 for details. The presented calculations 



make use of Newtonian gravity which leads to less compact 
neutron stars than General Relativity. This could have a 
moderate effect on the amount of ejecta mass and, since 
the /3-equilibrium Y e is determined by the density inside the 
initial neutron star, the electron fraction. As we will show 
below, however, the heavy nucleosynthesis is only affected 
once Y e reaches ~ 0.2 and the ejecta masses are (within the 
existing uncertainties) still consistent with compact binary 
mergers being a major r-process source even if the amount 
of ejecta is reduced by a factor of a few. Therefore, we ex- 
pect our major results to be robust. Nevertheless this point 
needs improvement in future work. 



2.3 Ejecta properties 

The merger dynamics is rather sensitive to asymmetries in 
the component masses, see Fig. [l] (for a more complete 
overview over the remnant structures see Fig. 1 in |Piran| 



et al. (2012)). The asymmetries enhance the degree to which 



the secondary is disrupted and they lead to larger ejecta 
masses and velocities. As already pointed out by Oechslin 
|et al.| ( |2007[ ), the ejecta stem from two different regions: 
a hot component from the interaction region between the 
stars (subsequently called "interaction-" or i-component for 
short) and a colder one where matter is simply flung out 
by tidal torques, see the first two panels in Fig. [2] We re- 
fer to this matter as "tidal" or t-component. As the mass 
ratio decreases, the primary only ejects small amounts of 
the i-component, while the fraction and overall amount in 
the t-component from the secondary increases. This differs 
from the findings of |Goriely et et\~] ( |2011| ), where even in 
unequal-mass cases most of the mass is lost from the pri- 
mary as i-component. For now we can only speculate about 
the origin of this difference: it could either be due the differ- 
ent treatment of self-gravity (conformal flatness vs Newto- 
nian) or, since the i-component comes from a region of very 
strong shear, it could also to some extent be impacted by 
the specifics of numerical method, such as implementation 
of articifical viscosity. 

As can be seen from the volume rendering in panel 3 of 
Fig. [2] the ejecta are extremely neutron-rich. In Fig. [3j left 
panel, we bin the mass fractions of the ejecta of all sim- 
ulations according to Y e . The values are clustered around 
Y e = 0.03 and dominated by t-component material. Since at 
high densities the V e -changing reactions are Pauli-blocked 
and the subsequent rapid expansion leads to a very fast tem- 
perature drop, the t-component had essentially no time to 
change its electron fraction. Therefore, it still carries the 
original Ye-value set by cold, high-density /3-equilibrium in- 
side the neutron star. The Y e of this matter therefore is a 
unique function of density which shows up as well-defined 
curve in the inset of Fig. [3] The scattered component in this 
inset is a combination of matter from the i-region and matter 
from the high Y e neutron star crust which is only marginally 
resolved in our simulations. The t-component of the ejecta 
not only starts within a narrow range of thermodynamical 
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Table 1. Overview over the performed simulations, the superscript + indicates that the primary is a black 
hole. 



Run 


mi [M ] 


m 2 [M ] 


7V S PH [10 6 ] 


tend H 




[M ] 


comment 


1 


1.0 


1.0 


1.0 


15.3 


7.64 x 


10 -3 




2 


1.2 


1.0 


1.0 


15.3 


2.50 x 


10 — 2 




3 


1.4 


1.0 


1.0 


16.5 


2.91 x 


10 -2 




1 


1.6 


1.0 


1.0 


31.3 


3.06 x 


10~ 2 




5 


1.8 


1.0 


1.0 


30.4 


> 1.64 x 


IO" 2 


secondary still orbiting 


5 


2.0 


1.0 


0.6 


18.8 


> 2.39 x 


10 — 2 


o t_ 1 1 cLl y oLlll \Jl L/lLlllti 


7 


1.2 


1.2 


1.0 


15.4 


1.68 x 


10~ 2 




8 


1.1 


1.2 


1.0 


13.9 


2.12 x 


10" 2 




g 


1.(5 


1.2 


1.0 


14.8 


3.33 x 


10 — 2 




10 


1.8 


1.2 


1.0 


21.4 


3.44 x 


10~ 2 




11 


2.0 


1.2 


0.6 


15.1 


> 2.95 x 


10~ 2 


secondary still orbiting 


12 


1.1 


1.4 


1.0 


13.4 


1.28 x 


IO" 2 




13 


1 . 6 


1.4 


1.0 


12.2 


2.36 x 


10 -2 




14 


1.8 


1.4 


1.0 


13.1 


3.84 x 


10~ 2 




15 


2.0 


1.4 


0.6 


15.0 


3.89 x 


IO" 2 




16 


1.6 


1.6 


1.0 


13.2 


1.97 x 


lO" 2 




17 


1.8 


1.6 


1.0 


13.0 


1.67 x 


10~ 2 




18 


2.0 


1.6 


0.6 


12.4 


3.79 x 


io- 2 




19 


1.8 


1.8 


1.0 


14.0 


1.50 x 


lO" 2 




20 


2.0 


1.8 


0.6 


11.0 


1.99 x 


10~ 2 




21 


2.0 


2.0 


0.2 


21.4 


1.15 x 


10" 2 




22 


5.0 + 


1.4 


0.2 


138.7 


2.38 x 


10~ 2 


nsbh 


23 


10.0 + 


1.4 


0.2 


139.0 


4.93 x 


io- 2 


nsbh 



parameters, Y e =0.04 ±0.02, p= (1.4 ±0.5) x 10 1 ' gem 
corresponding to the outer ~ 2 km of the original neutron 



star (see, for example, Fig. 1 in Rosswog et al. (2012 1), but 
also continues to expand adiabatically and without shocks. 
Therefore, all fluid elements in tidal ejecta share approxi- 
mately the same thermodynamic history and provide robust 
conditions for a uniform nucleosynthesis. 
The amounts of ejected matter range from 7.6 x 10" 3 to 
3.9 x 10" 2 M Q (4.9 x 10" 2 M in the nsbh case). We see a 
clear trend for decreasing mass ratios to eject more matter, 
but for the equal-mass systems no unique trend is found with 
respect to the total binary mass. For at least small ns mass 
asymmetries, i.e. for the most common cases, we provide 
a fit formula based on the dimensionless mass asymmetry 
parameter 



4mim2/(mi ± 7712) 



(1) 



extreme mass ratios^ The fit formula 
m e j(mi,m2) = (rn-i + 1x12) I A — Brj — 



C 



l±7/ 3 /0- S 



(2) 



with A = 0.0125, B = 0.015, C = 0.0083 and a = 0.0056, 
shown as solid line in Fig. [3] right panel, provides a rea- 
sonably good approximation to the simulation results. For 
unequal-mass binaries m e j/(mi +m,2) tends to a minimum 
when approaching r\ — > and a maximum of w 0.12 is 
reached near 77 = 0.02. 



3 NUCLEOSYNTHESIS 

3.1 Hydrodynamic trajectories 

Our nucleosynthesis calculations are performed with, a large 
reaction network ( |Winteler|2012| |Winteler et al.||2012[) t hat 
is based on the BasNet network ( Thielemann et al.|2011 1. It 
includes over 5800 isotopes from nucleons up to Z = 111 be- 
tween the neutron drip line and stability. The reaction rates 



are from the compilation of Rauscher & Thielemann ( 2000 1 



Contrary to the mass ratio q, rj is symmetric with respect 
to both masses and varies only in a finite range from to 
1 (with for a equal-mass system and rj approaching 1 for 



2 A similar parameter rj = rn\m2lirn\ + m.2) 2 was used in de- 
riving an empirical formula for recoil velocities of coalescing non- 
spinning black holes (iGonzalez et al. 120071 |Fitchett|1983| 
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Figure 2. The ejecta come from two different regions (left and middle): a hot interaction region between the stars where matter is 
ejected by hydrodynamic effects and a colder region that is flung out by tidal torques. The fraction of the latter material increases with 
the asymmetry in the stellar masses (left f.4 + 1.4 Mq, right 1.8 + 1.4 Mq). Right: volume rendering of the Y e distribution (1.4 + 1.2 
Mq at t= 8.09 ms), only matter below the orbital plane is shown. 
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Figure 3. Left panel: Distribution of Y e in the selected subset of ejected particles in all ns 2 and nsbh simulations, binned by mass. 
The inset shows how Y e is correlated with the initial density in each particle before the merger. Majority of points on the Y e -p diagram 
trace the neutron star composition in the layer below the crust, where the Y e (p) reaches minimum. As a consequence, Y e in the ejecta is 
narrowly distributed around Y e ~ 0.04. Right panel: Fit of the ejected mass normalized by the total mass of the binary as a function of 
the asymmetry parameter rj (see main text for definition). Coefficients of the fit are: A = 0.0125, B = 0.015, C = 0.0083, and a = 0.0056. 
The arrows indicate the lower limits on the ejected mass for the three simulations in which the secondary is still not fully disrupted. 



for the finite range droplet model (FRDM) (Moller et al. 
1995 ) and the weak interaction rates (electron/positron cap- 



tures and /3-decays) are the same as in Arcones & Martinez- 



Pinedo (20111. In addition, neutron capture and neutron- 



induced fission rates of Panov et al. (20101 and /3-delayed 



fission probabilities as described in Panov et al. (2005) are 
used. 



The thermodynamical conditions are taken from hydrody- 
namical trajectories of individual SPH particles. We calcu- 
late nucleosynthetic yields for 30 representative trajectories 
in each of the 21 ns 2 merger simulations, and for 20 trajecto- 
ries in the two nsbh merger simulations. These trajectories 
are a fair representation of the total ejected mass, see Tab.[l] 
Once matter has been ejected from the central high-density 
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Figure 4. Top: density and temperature evolution for a bundle 
of trajectories, color-coded by density at T = 10 GK. Middle: 
resulting final abundances distribution. Their averaged distribu- 
tion is shown in a black dashed line, and a bold red line represents 
abundances for a trajectory without heating. All trajectories rep- 
resent a subset from the standard 1.4-1.4 Mq merger. Bottom: 
distribution of abundances for a variety of different (ns 2 and nsbh) 
merger cases. All different astrophysical systems yield essentially 
identical resulting abundances. 
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Figure 5. Evolution of the neutron density N n (top), aver- 
age neutron separation energy {S n ) (middle) and average pro- 
ton number (Z) (bottom) for five trajectories from run 12 with 
different initial densities (in the range 10 12 — 10 13 g/cm 3 ). 



region, both temperature and density decrease rapidly while 
the electron fraction Y e can, in principle, start to increase 
due to /3-decays and positron captures. A small fraction of 
the trajectories (less than 0.008, see Fig. [3] left panel) suf- 
fers shocks during their later evolution during which their 
Y e increases up to ~ 0.25, but for the bulk of trajectories the 
Ve-values remain essentially constant. Some of the higher Im- 
material also comes from the (marginally resolved) neutron 
star crust material. Although heavily dominated by values 
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around 0.03, our ejecta electron fractions show a broader 



range than those of Goriely et al. (20111, probably due to 



In our calculations the low V e results in an initial composi- 
tion consisting of neutron-rich nuclei (Z > 20, on the neu- 



our higher resolution (crust material) and their omission of 
weak interactions. 

We start our nucleosynthesis calculations when the temper- 
ature of the trajectory has dropped below To = 10 GK. 
The initial composition is given by nuclear statistical equi- 
librium and it is further evolved with the complete network. 
Because most of the trajectories from hydrodynamical sim- 
ulations terminate after t^n ~ 10 — 20 ms, we extrapolate 
according to a free expansion for the density and an adia- 
batic expansion law for the temperature: 



P{t) = Pfln 



T(t)=T[S an ,p(t),Y e (t)}. 



where pa n and Sfi n are the density and entropy at ta n , while 
the temperature is calculated at each point of time from the 
equation of state ( Timmes & Swcsty 2000 1 . 
Contrary to neutrino-driven winds where the ejecta are dom- 
inated by alpha particles and only few seed nuclei partici- 
pate in the r-process, essentially all matter in neutron star 
mergers undergoes r-process. Therefore the r-process energy 
generation is non-negligible and leads to a substantial tem- 
perature increase (see e.g., Freiburghaus et al. 1999||Metzger| 
et al. 2010 | Goriely et al. 2011 1 that could affect the late-time 
evolution ( |Metzger et al.||2010 |. We calculate the r-process 
heating in a post processing step and consider its influence 
on the nucleosynthesis by modifying the temperature. To 
account for neutrino energy losses associated with /3-decays, 
we introduce a heating efficiency parameter e t h which mea- 
sures the fraction of nuclear power which is retained in the 
matter. Metzger et al. (20101 argue that this fraction must 
be e t h ~ 0.25 — 1. As a default, we use e t u — 0.5, but below 
we explore how (in)sensitive the results are to this choice. 
For a given trajectory we take its initial entropy Sh y d(i) 
and increment it by tthSnuc(t) , following the prescription of 



Freiburghaus et al. ( 1999 1 



3.2 Final abundances 

The r-process nucleosynthesis in compact binary mergers is 
extremely robust, all 23 cases of Tab. [I] deliver essentially 
identical abundance distributions. We illustrate the origin 
of this robustness using the reference case of two neutron 
stars with 1.4 M each (run 12 in Table[T|. Fig.|I]shows the 
density and temperature evolution for all trajectories (upper 
panel) and their final abundances. Between different trajec- 
tories of a single run there is only a very small spread in 
the resulting abundances, although the evolution of temper- 
ature is not always the same. The resulting final abundances 
(sum over all ejecta trajectories for each run) are practically 
the same for all cases. We show a representative selection of 
the results from different runs in the bottom panel of Fig. [4] 
The robustness of the abundances is mainly due to two fac- 
tors: 1) the r-process path always reaches the drip line be- 
cause of the extremely neutron-rich conditions and 2) there 
are several fission cycles (similar to what was observed in 
Goriely et al.|[20TT I . These conditions lead to the evolution 
of the neutron density shown in Fig. [5] upper panel. The 
huge values of the neutron density are due to the (narrow 
range of) very low Y" e -values. 



tron drip line) and neutrons. Differences with Goriely et al. 
(2011) in lighter heavy elements are not particularly relevant 
because these can also be produced in other astrophysical 
events, and the overall yields of the elements with A < 120 
that we obtain are much smaller compared to the heavy ro- 



bust r-process elements (similar to Goriely et al. (20111) 



Note that the contribution of the high-Ke trajectories to the 
final abundances is negligibly small. 

In order to understand the evolution of the r-process path 
we monitor the average neutron separation energy (middle 
panel Fig. [SJ defined as 

, q , J2z.ASn(Z,A)Y(Z,A) 

{n) ~ H Z , A Y{Z,A) ' (6) 

with Sn(Z, A) and Y(Z,A) being the neutron separation 
energy and abundance of the nucleus (Z, A) . The average 
separation energy decreases when matter moves away from 
stability and it is, by definition, zero at the neutron drip 
line. The second panel of Fig. [5] shows that the average 
separation energy is initially below ~ 1 MeV, which indi- 
cates that the r-process path proceeds along the neutron 
drip line. The average proton number increases to Z — 40 
at t ~ 10~ 2 s where the neutron separation energy reaches 
a maximum. This local maximum occurs when the magic 
number TV = 82 is overcome. Whenever the r-process path 
reaches a neutron magic number, it moves closer to the line 
of ^-stability (i.e., larger S n values) by increasing Z with- 
out changing TV. After the matter flow passes TV = 82 (here 
around Z = 40), the (S n ) decreases because the path gets 
again further away from the /3-stable region. This contin- 
ues until the next magic number, TV = 126, is encountered 
(corresponding to the minimum of (S n ) between 0.01s and 
0.1s). Around t = 0.1s and Z = 60, TV = 126 is also over- 
come as indicated by the second maximum of (S„). After 
this point the oscillations in (S n ) are due mainly to fission 
as can be seen by the behavior of (Z). Note that FRDM 
predicts a magic number at TV = 184, consequently the path 
reaching this point may lead to oscillations in (S n ). There- 
fore, the quantity (Z) is a nuclear fission indicator that is 
better suited than the average mass number ( Goriely et al. 
|2011[ ), since it can decrease only through fission reactions, 
while the average mass number can also decrease due to 
photo-dissociations. The maximum (Z) corresponds to the 
moment when a significant amount of matter has reached 
the region where fission becomes important. The daughter 
nuclei resulting from fission capture neutrons move first to- 
wards the drip line and then to higher Z where fission acts 
again. In that way several fission cycles occur and lead to 
oscillations in (Z). In Fig. [5] one can distinguish at least 
three fission cycles. The final increase of (S n ) is due to the 
r-process freeze-out at t w Is when the neutron-to-seed ratio 
drops below unity and the matter /3-decays to stability. We 
have used our reference case of two 1.4 Mq neutron stars, 
run 12, for illustration purposes, but behaviour is very sim- 
ilar in all of the other cases. 

The major result of this study is that the r-process abun- 
dances for compact binary mergers do not depend on the 
astrophysical parameters of the merging system. They do 
depend, however, on the not-so- well-known properties of nu- 
clei close to the neutron drip line. The most critical inputs 
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Figure 6. Distribution of final abundances in the vicinity of the second and third r-process peaks. Left: for the case without heating 
and for three different values of an heating efficiency e^; Right: for two types of fission fragments mass distribution. Overplotted black 
dots represent solar abundances. 
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Figure 7. Left: postprocessed temperature as a function of pressure for the cases without nuclear heating (dashed line) and with 
heating efficiencies = 0.1, 0.5 and 0.9. Right: nuclear heating power as a function of time, along with a power-law fit <x t~ a with 
index a = 1.3. The small black triangles on each trajectory are neutron freeze-out points (the points in which Y n /Y scc ^ = 1). 



are the nuclear masses, which determine the location of the 
drip line and the evolution of the path, and the fission bar- 
riers and yields distribution, which are responsible for the 
abundances between 120 < A < 170. We explicitely explore 
how sensitive our results are to the details of our nuclear 
physics input. In Fig. [6] left panel, we show the effect of 
including heating with different efficiencies of the profile of 
final abundances. In all cases the second and third r-process 
peak are robustly reproduced, but the regions behind the 
peaks are less populated if heating is included. In the middle 
panel we vary the fraction of nuclear energy that is retained 
in the matter, eth- As long as eth is above some threshold 
value (~ 1%) the resulting abundances are independent of 
the exact choice. In the right panel we explore the sensitiv- 
ity to the fission prescription. We once follow [Panov et al.| 
( 2001 1 assuming only two fission daughter nuclei and, alter- 



natively, we adopt the fission treatment and yield distribu- 



tion of Kodama & Takahashi ( 1975 I. When using the latter, 



the trough after the second peak completely fills up while 
for the Panov prescription the fission products end up en- 
tirely in the second peak. Note that the used prescription 
of Panov uses only two daughter nuclei while he also pro- 
vides yield distributions in a more recent paper ( Panov et al. 
|2008[ ). Therefore the deviation of our results from the solar 
system abundances are sensitive to the physics of extremely 
neutron-rich nuclei. 

In Fig. [7] left panel, we show the evolution of a trajectory 
in the pressure-temperature plane for different heating effi- 
ciencies, the right panel shows the generated nuclear power. 
The early evolution is very similar for all four cases because 
the internal energy is high and the contribution from the 
generated heat is negligible. The latter becomes the domi- 
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nant part of the internal energy after about 20 ms when the 
temperature and density decrease enough for r-process to 

dl999b we find that the 



begin. Similar to Freiburghaus et al.| 
generated nuclear energy raises the temperature by almost 
an order of magnitude compared to the case without heat- 
ing. Matter essentially stays in (n, 7)-(7, n) equilibrium until 
the neutron-to-seed ratio drops below unity ("freeze-out"). 
Note that all non-zero heating efficiencies produce trajecto- 
ries that lie on the P — 1 /3aT 4 -line at freeze-out and there- 
fore follow a similar behaviour in their decay to stability. The 
reason is that for our extremely neutron-rich conditions the 
main contribution to the pressure comes from neutrons and 
photons, while the electron pressure is orders of magnitude 
smaller. Therefore, at the freeze-out point when most of the 
neutrons are captured the system is radiation pressure dom- 
inatec0 

Later the system evolution enters the radiation-dominated 



expansion which was studied in Li & Paczyhski (19981. In 



this regime the generated nuclear heat comes mainly from 
/3-decays and can overall be well approximated by a power 
law, while before the freeze-out the generated heat stays ap- 
proximately constant (see Fig. We suggest the following 
fit which smoothly interpolates between a constant value 
and a power law: 



e(t) = to 



1 t - 

— arctan 

7r a 



to 



0J> 



(4) 



with e = 2 x 10 18 erg/(gs), to = 1.3 s, a = 0.11 s, and 
a = 1.3. This formula yields a good fit to all heating histo- 
ries. Similarly to the trajectory considered in Mctzgcr ct al. 
(20101, in our case starting from ~ 10 s the nuclear heating 



power becomes essentially independent from the initial con- 
ditions. The deviations from a power law that are visible in 
the right panel of Fig.JTJat t ~ 10 _1 d and t ~ 20 d are due 
to the decay of individual elements. 

To explore the dependence of the final abundances on the 
electron fraction Y e , we show in Fig. [8] the resulting abun- 
dances for some of the few trajectories with higher Y e -values 
(colour; for comparison the results from the dominating low- 
Ye trajectories are shown in gray). For Y e up to ~ 0.15 the 
results hardly deviate from each other, only at and above 
this value do the abundances below A ~ 110 start to devi- 
ate. Substantial deviations for the abundance pattern of the 
heaviest nuclei with A > 110 only occur for Y e > 0.2, since 
here the Y e increase is triggered by shocks at low density, 
the neutron density drops earlier and the neutron captures 
freeze out. 



4 SUMMARY AND DISCUSSION 

In this study we have re-examined the question of heavy 
element nucleosynthesis in compact binary mergers. Our 
study adds several new aspects in comparison with earlier 
work. First, we systematically cover the plausible ns 2 
parameter space with masses from 1.0 to 2.0 M© in 21 
simulations. Despite the long history of ns 2 simulations we 

3 At higher temperatures the radiation pressure must necessar- 
ily include contribution from electron-positron pairs l |Witti et al.| 




1994| |Farouqi et al.|2uT0 l 
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Figure 8. Dependence of the final abundances on the electron 
fraction. Shown are the final abundances of some trajectories with 
high Y e (in colour), they are compared to the low-Ys trajectories 
(gray) that dominate the ejecta. Below Y e = 0.1 two or more 
fission cycles occur, one or two are realized for 0.1 < Y e < 0.22 
and none for Y e > 0.22. Since high-Y e trajectories only repre- 
sent a very small mass fraction, their contribution to the overall 
abundances is minuscule. 



are not aware of any study that has systematically explored 
such a wide range of ns masses. We complement these cases 
with mergers of a 1.4 M neutron star with black holes 
of 5 and 10 M©. Second, since we use a nuclear equation 
of state with neutron stars in initial cold /3-equilibrium 
and include Y e -changing weak interactions, we overcome 
the Ye-ambiguity that has plagued all previous studies. We 
subsequently find the final abundances by following the nu- 
cleosynthesis along a large set of hydrodynamic trajectories 
with a state-of-the-art nuclear reaction network. 

Consistent with earlier studies, we find a very robust 
r-process and final abundances that are in good agreement 
with the heavy solar system abundance pattern. The major 
new result is that the final pattern is extremely robust 
accross the whole parameter space: all 21 ns 2 -merger 
and the two nsbh-merger cases yield practically identical 
nucleosynthesis outcomes. The major reason for this unique 
abundance pattern is the extreme neutron richness of the 
ejecta, (Y e ) w 0.04. Consequently, in each case the r-process 
path meanders along the neutron drip line and matter 
undergoes several fission cycles, so that the abundances are 
determined entirely by nuclear rather than by astrophysical 
properties. As a corollary, the poorly known nuclear prop- 
erties near the neutron dripline do have an impact on the 
resulting abundance pattern. We find some dependence on 
the used mass formula and on the distribution of the nuclear 
fragments after fissioning. Nevertheless, the second and 
third r-process peaks are robustly reproduced, without any 
"tuning" of the nuclear physics input the overall agreement 
with the solar system r-process pattern is good. R-process 
matter lighter than ~ A ~ 120, however, is substantially 
under-produced with respect to the solar system pattern. 
The few high-Y e trajectories produce different abundance 
patterns, see Fig. [SJ but there is too little mass in this 
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Figure 9. Relevance of the ejected mass for the r-process en- 
richment of the galaxy. The mass asymmetry parameter rj, see 
Eq. is used here to label the different simulations. Adopt- 
ing the 95% confidence intervall for the galactic double neutron 
star merger rate from | |Kalogera et al.|2004| |, we calculate which 
galactic enrichment rates would result for the ejected matter of 
our simulations (lines with squares at ends; our reference case of 
two 1.4 Mq neutron stars is shown in magenta). For comparison 
we also plot the estimated galactic enrichment rate | |Qian|2000| |. 
All the ejected masses are consistent with being a major con- 
tributor to galactic r-process, even if the total amount per event 
should be decreased by a factor of a few (say, due to equation of 
state or GR effects). 



material to leave a noticable impact in the resulting average 
distribution. 

One may wonder whether captures of neutrinos from the 
central remnant, which are not included in our simulations, 
could change the Y e of the ejecta noticeably. For two 
reasons, we do not think that this is the case. First, the 
ejecta are launched practically at first contact, see Fig. [2j 
with (v) > 0.1c, ( Piran et al.|2012 l. The neutrino emission, 
however, only becomes relevant after a disk has formed 
which takes ~ 7 more milliseconds (e.g. Fig. 7 in Rosswog 



et al. (2012)) at which time the solid angle under which the 
escaping matter is seen has decreased already by a large 
factor. Second, the ejecta are concentrated in the orbital 
plane while the neutrinos are emitted preferentially perpen- 



dicular to it, see Fig. 12 in Rosswog & Liebendorfer (20031 



and Fig. 9 in Dessart et al. (20091. We therefore expect 



Ye-changes due to neutrino captures to be negligibk^J 
The dynamic ejecta masses are large enough to make a 
substantial contribution to the chemical enrichment of the 
galaxy with r-process material. Adopting the observation- 
based double neutron star merger rates of |Kalogcra et al.| 
( 2004b , -Rons = SS.Otf^ 1 yr" 1 , we calculate by which 



average rate the galaxy would be r-process enriched from 
neutron star mergers. Although the overall ejecta masses 
vary by more than a factor of five between different cases, all 

4 This may be different for likely emerging neutrino-driven winds 
that might also contribute to the nucleosynthesis of compact bi- 
nary merger events. 



of them are perfectly consistent with the estimated galactic 
enrichment rate ( Qian|20"00 l, see Fig. [9] The reference case 
of 2 x 1.4 Mq neutron stars (shown in magenta), nicely 
brackets this value. This also implies that the relevance of 
compact binary mergers for r-process remains unaffected 
even if the typical ejecta amount per event should be 
decreased by a factor of a few, say, due to equation of state 
or GR effects. 

Our findings of ejecta masses in the required range, 
extremely robust abundance patterns for the heaviest 
r-process elements, but underproduction of the lighter 
r-process elements, make compact binary mergers the 
major candidate for the observed, unique, heavy r-process 
component (Sneden et al. 20081. At the same time this 



means, that at least one, but possibly more additional 
r-pocess sites must contribute to the cosmic chemical 
evolution of lighter r-process elements. 

An important question is how early the r-process yields 
from compact binary mergers would be available for cosmic 
enrichment. The inspiral time from the birth of a compact 
binary system to a coalescence is very sensitive to both 



the initial sparation and the orbital eccentricity (Peters 



& Mathews 1963 19641. The initial separation is set by 



binary evolution processes which are beyond the scope of 
this paper, but it is worth re-iterating that a binary system 
that survives the two SN explosions that are required to 
form the system in the first place will in most cases possess 
a large orbital eccentricity. As illustrated in Fig. [10] short 
lifetimes can be achieved with plausible orbital parameters. 
For example, an initial semi-major axis of 1 R and an 
eccentricity of 0.9 lead to a merger only 1 Myr after binary 
formation (for comparison we note that the projected 
semi-major axis of the double pulsar PSR J0737-3039 



A+B (Burgay et al. 2003 Kramer & Stairs 20081 is 0.6 



Rq). Some binary population synthesis studies (Dewi & 
[Pois][2003l |Belczynski et al.||200u] ) find indeed short-lived 
channels that coalesce essentially at birth, so an early 
enrichment from compact mergers is at least plausible. 
Earlier studies ( Argast et al.|2004| disfavoured neutron star 
mergers as dominant galactic r-process source. Based on 
our hydrodynamic plus nucleosynthesis studies, however, 
we have shown that compact binary mergers are excellent 
production sites for the unique, heavy r-process component. 
We leave the question of whether or to which extent they are 
consistent with cosmic chemical evolution for future studies. 
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